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1.  INTRODUCTION 


The  usual  mixed  linear  model  discussed  in  Che  literature  on 
variance  components  is 

Y - X0  + U,*,+...+U  * + c (1.1) 

li  p p 

where  X,  VI,,..., U are  known  matrices,  0 is  a fixed  unknown 

1 P 

vector  parameter  and  4^,...,$^,  c are  unobservable  random 
variables  (r.v. ’s)  such  that 


E(c)  - 0.  E<*t)  - 0,  E^j)  - 0,  i*J,  E(e*p  - 0 


«**•>  - <&„•  e<v;>  - "k,- 


(1.2) 


2 2 2 

The  unknown  parameters  oq,  o^,...,  o are  called  variance  components. 


Some  of  the  early  uses  of  such  models  are  due  to  Yates  and 
Zacopancy  (1935)  and  Cochran  (1939)  in  survey  sampling,  Yates  (1940) 
and  Rao  (1947,  1956)  in  combining  intra  and  interblock  information 
in  design  of  experiments,  Fairfield  Smith  (1936),  Henderson  (1950), 
Panse  (1946)  and  Rao  (1953)  in  the  construction  of  selection  indices 
in  genetics,  and  Brownlee  (1953)  in  industrial  applications.  A 
systematic  study  of  the  estimation  of  variance  components  was  under- 
taken by  Henderson  (1953)  who  proposed  three  methods  of  estimation. 

The  general  approach  in  all  these  papers  was  to  obtain  p+1 
quadratic  functions  of  Y,  say  Y'Q.Y,  1 ■ l,...,p+l,  which  are 
invariant  for  translation  of  Y by  Xa  where  a is  arbitrary, 
and  solve  the  equations 


1.1 


(1.3) 


Y’Q.Y  - E(Y'Q.Y)  - a,  o2  + a . to2+. . .+a4  o2 
l 1 lo  o il  1 ip  p 
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The  method  of  choosing  the  quadratic  forms  was  intuitive  in  nature 
(see  Henderson,  1953)  and  did  not  depend  on  any  stated  criteria  of 
estimation.  The  entries  in  the  ANOVA  table  giving  the  sums  of 
squares  due  to  different  effects  were  considered  as  good  choices 
of  the  quadratic  forms  in  general.  The  ANOVA  technique  provides 
good  estimators  in  what  are  called  balanced  designs  (see  R.  L.  Anderson, 
1975  and  R.  L.  Anderson  and  P.  P.  Crump,  1967)  but,  as  shown  by 
Seely  (1975)  such  estimators  may  be  inefficient  in  more  general  linear 
models.  For  a general  discussion  of  Henderson's  methods  and  their 
advantages  (computational  simplicity)  and  limitations  (lack  of 
uniqueness,  inapplicability  and  inefficiency  in  special  cases)  the 
reader  is  referred  to  papers  by  Searle  (1968,  1971),  Seely  (1975), 

Olsen,  Scaly  and  Birkes  (1976)  and  Harville  (1977,  p.  335). 

A completely  different  approach  is  the  ML  (maximum  likelihood) 

method  Initiated  by  Hartley  and  Rao  (1967).  They  considered  the 

likelihood  of  the  unknown  parameters  8,  o",...,o  based  on  observed 

o p 

Y and  obtained  the  likelihood  equations  by  computing  the  derivatives 

of  likelihood  with  respect  to  the  parameters.  Patterson  and  Thompson 

(1975)  considered  the  marginal  likelihood  based  on  the  maximal 

invariant  of  Y,  i.e.,  only  on  B'Y  where  B is  a matrix  orthogonal 

to  X and  obtained  what  are  called  restricted  maximum  likelihood  (RML) 

equations.  Harville  (1977)  has  given  a review  of  the  ML  and  RML 
methods  and  the  computational  algorithms  associated  with  them. 

0 

One  criticism  of  the  Ml.  estimators  is  that  they  may  be  heavily 
biased  so  that  some  caution  is  needed  when  they  are  used  as  estimates 


i t 


i im  n> mw* 


. 


of  Individual  parameters  for  taking  decisions  or  for  using  them 
in  tha  placa  of  true  values  to  obtain  an  efficient  estimate  of  B. 

Tha  problem  is  not  accute  if  the  exact  distribution  of  the  ML 
estimators  is  known,  since  in  that  case  appropriate  adjustments 
can  be  made  in  the  Individual  estimators  before  using  them.  The 
general  large  sample  properties  associated  with  ML  estimators  are 
misleading  in  tha  absence  of  studies  on  the  orders  of  sample  sizes 
for  which  these  properties  hold  in  particular  cases. 

The  bias  in  RML  estimators  may  not  be  large  even  in  small 
samples.  But  both  ML  and  RML  estimators  are  functions  of  B'Y, 
the  maximal  invariant  of  Y,  and  there  are  important  practical  cases 
where  reduction  of  Y to  B'Y  results  in  non-identif lability  of 
individual  parameters,  in  which  case  neither  ML  nor  RML  is  applicable. 

The  details  ara  given  in  Section  5. 

Rao  (1970,  1971a,  1971b,  1972,  1973)  proposed  a general  method 
called  MINQE  (Minimum  norm  quadratic  estimation)  the  scope  of  which 
has  been  extended  to  cover  a variety  of  situations  by  Focke  and  Dewess 
(1972),  Kief fa  (1975,  1976,  1977a, b,  1978,  1979),  J.  N.  K.  Rao  (1973), 

Fuller  and  Rao  (1978) , fbduri  Rao  and  Chaubey  (1978),  Pukclsheim  (1977, 

1978),  Slnha  and  Wleand  (1977)  and  Rao  (1979).  The  method  is  applicable 
to  a general  linear  model 

Y - X8  + c,  E(ce')  - 0^+. . ,+6pVp  (1.4) 

where  no  structure  need  be  Imposed  on  c and  no  restrictions  are 
placed  on  0^  or  V^.  (In  the  model  (1.1),  0^  >_  0 and  are 

non-negative  definite). 

* In  the  MINQE  theory,  we  define  what  is  called  a natural  estimator 


1.3 


of  a linear  function  f'6  of  6 In  term*  of  tha  unobaervable  r.v. 
e In  (1.4),  aay  t'Ne.  Then  tha  estimator  Y'AY  In  terms  of  tha 
observable  r.v.  Y Is  obtained  by  minimising  the  norm  of  the  difference 
between  the  quadratic  forms  e'Ne  and  Y'AY  • (X0  + c)'A(XB  + t). 

The  universality  of  the  MTNQfmethod  arises  from  the  following 
observations: 

(a)  It  offers  a wide  scope  in  the  choice  of  the  norm  depending 
on  the  nature  of  the  model  and  prior  Information  available. 

(b)  One  or  more  restrictions  such  as  Invariance,  unbiasedness 
and  non-negative  definiteness  can  be  placed  on  Y'AY 
depending  on  the  desired  properties  of  the  estimators. 

(c)  The  method  is  applicable  is  situations  where  ML  and  RML  fail. 

(d)  There  is  an  automatic  provision  for  Incorporating  available 
prior  information  on  the  unknown  parameters  6 and  6. 

(e)  Further,  ML  and  RML  estimators  can  be  exhibited  as  Iterated 
versions  of  suitably  chosen  MINQts. 

(f)  The  MINQC equation  provides  a natural  numerical  algorithm  for 
computing  the  ML  or  RML  estimator. 

(g)  For  a suitable  ch->i  ro  of  the  norm,  the  MINQ  estimators  provide 
minimum  variance  estimators  of  0 when  Y is  normally 
distributed. 

It  has  been  mentioned  by  some  reviewers  of  the  MINQE  theory  that 
the  computations  needed  for  obtaining  the  MINQ  estimators  are  somewhat 
heavy.  It  Is  true  that  the  closed  form  expressions  given  for  MINQE' s 
contain  Inverses  of  large  order  matrices,  but  they  can  be  computed  In 
a simple  way  In  special  cases  that  arise  in  practice.  The  computations 
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In  such  cases  are  of  the  same  order  of  magnitude  as  obtaining  sums 
of  squares  in  the  ANOVA  table  appropriate  for  the  linear  model. 
Perhaps,  some  research  is  needed  in  developing  simple  numerical 
algorithms  in  more  complicated  cases.  It  is  certainly  not  true 
that  the  computation  of  MLE  or  RMLE  is  simpler  than  that  of  MINQE. 
Both  may  have  the  same  order  of  complexity  in  the  general  case. 


2.  MODELS  OF  VARIANCE  AND  COVARIANCE  COMPONENTS 


2.1  General  Model 

There  la  e large  variety  of  models  of  variance  end  covariance 
component a used  In  reaaarch  work  In  biological  and  behavioral  sclencea. 
They  can  all  be  considered  In  a unified  frame  work  under  a general 
Gauea-Markoff  (CM)  modal 


Y - XB+  € (2.1.1) 

where  Y la  n-vector  random  variable,  X la  n x m matrix,  g is 
m-vector  parameter  and  c la  n-vector  variable.  The  models  differ 
mainly  In  the  structure  imposed  on  c.  The  most  general  formulation  Is 


E(e)  - 0 (2.1.2) 

d(c)  - e.v,+...+ev  - v(e>  - va  (2.1.3) 

ii  p p e 

where  D atanda  for  the  dispersion  (variance  covariance)  matrix, 

A'  » (0^,...,0p)  is  unknown  vector  parameter  and  V^,...,V^  are 


known  symmetric  matrices.  We  let  0cR°  and  0cF  (open  set)  c RP 
such  that  V(0)  >_  0 (l.e.,  nonnegative  definite).  In  the  representation 

(2.1.3)  wa  have  not  Imposed  any  restriction  such  as  0^  > 0 or  is 

nonnagatlve  definite. 

It  may  be  noted  that  any  arbitrary  n x n dispersion  matrix 
9 ■ (S^j)  can  be  written  In  tha  form  (2.1.3) 

rr6iJVlJ  (2.1.4) 

Involving  a maximum  of  p ■ n(n-l)/2  unknown  parameters  0^  and 
known  matrices  V^ , but  in  model*  of  practical  interest  p has  a 


2.1 


relatively  email  value  compared  to  n. 


2.2  Variance  Components 

A special  case  of  the  variance  components  model  is  when  e has 
the  structure 


c - UJ.+...+U* 

11  p p 

where 

Is  n x m^  given  matrix  and 

that 

E^)  - 0.  E(4t*p  - 0 lf«J 

(2.2.1) 

1 iB  m^-vector  r.v.  such 


In  such  a case 


V (0)  - e1vi+-*‘+epVp  (2.2.3) 

2 

where  ® *nd  Moat  of  the  models  discussed 

In  literature  are  of  the  type  (2.2.1)  leading  to  (2.2.3). 

In  the  classical  regression  model  p ■ 1.  Other  examples  are  one 
and  two  way  classification  models  with  fixed  and  random  effects. 

The  complete  GM  model  when  c has  the  structure  (2.2.1)  is 

Y - X0  + U141+...+Up4p 

E^)  - 0;  E^fj)  - 0 1 *){  E(^)  - a2Lln  . (2.2.4) 

The  associated  statistical  problems  are: 


(a) 

Estimation  of 

6 

(b) 

Estimation  of 

2 

i 1 m 1» • • • »P 

(c) 

Estimation  of 

1 • 1» « • * »P» 

2.2 


k 


The  last  problem  arises  in  Che  construction  of  selection  Indices  in 
genetics,  end  some  early  papers  on  the  subject  providing  a satisfactory 
solution  are  due  to  Fairfield  Smith  (1936),  anse  (19A6)  based  on  an 
idee  suggested  by  Fisher,  and  Henderson  (1950).  A theoretical 
Justification  of  the  method  employed  by  these  authors  and  associated 
tests  of  significance  are  given  in  Rao  (1953). 

A particular  case  of  the  model  (2. 2. A)  is  where  it  can  be  broken 
down  into  a number  of  submodels 

Y,  - X. 6 + c. , . . . ,Y  - X B + e (2.2.6) 

ill  p p p 

where  is  n^-vector  variable  and 

E(ct)  - 0,  E(Cle’)  - etIn  , E(Clcp  - 0.  (2.2.7) 

Mote  that  the  B parameters  are  the  same  in  all  submodels,  and  in 

some  situations  the  design  matrices  X. ,...,X  may  also  be  the  same. 

1 P 

The  model  (2.2.6)  with  the  covariance  structure  (2.2.7)  is  usually 
referred  to  as  one  with  "heteroscedastic  variances"  and  the  problem 
of  estimating  B as  that  of  estimating  a "common  mean". 

2.3  Variance  and  Covariance  Components 

We  assume  the  same  structure  (2.2.1)  for  c but  with  different 
covarlancaa  for  the  ^'s 

E^j)  - 0,  E(ft*')  - Ait  i - l,...,k 

E(*t4j)  - <JjIm  , i - k + l,...,p  (2.3.1) 

E(^i^J ) “ °*  1 * J 
2.3 


1 

H 


i 


>»■  -A.-. 
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leading  to 

V(8)  - U1A1u;+...+«k*kU-  ♦ (2.3.2) 

where  >_  0.  In  some  practical  problems  are  all  the  same  and 

2 

there  Is  only  one  o In  which  case  (2.3.1)  becomes 


V(0)  - U1AUj+...-HJkAU£  + o I.  (2.3.3) 


2.4  Random  Regression  Coefficients 

This  Is  a special  case  of  the  variance  and  covariance  components 
model  considered  In  2.3  where  c has  the  structure 

c - + *2.  E(^*p  - A,  E(*2*')  - o2I  (2.4.1) 

th(  compounding  matrix  for  ^ being  the  same  as  for  B leading  to 
the  CM  model 


y - xe  + x$L  + $2 


D(e)  - XAX'  + oz I. 


(2.4.2) 


In  general,  we  have  repeated  observations  on  the  model  (2.4.2) 
with  different  X’s 

T1  • V + Vll  + *21>  1 ■ 1 ‘ 


leading  to  the  model 


Y - XB  + € 


2.4 


(2.4.3) 


— - 


(2.4.4) 


with 


X]AX1  + o I 


0(c)  - 


(2.4.5) 


XtAXt  + 01 


all  the  off  diagonal  blocks  being  null  matrices.  A discussion  of 
such  models  is  contained  in  Fisk  (1967),  Rao  (1965,  1967),  Swamy  (1571), 
and  Spjotvoll  (1977). 


2.5  Intraclass  Correlation  Model 

We  shall  illustrate  an  intraclass  correlation  model  with  special 
reference  to  two  way  classified  data  udi  ~ cM 


»ijk- 1 ' 1 1 ‘ 1 k • 1 '• 


(2.5.1) 


We  write 


VWuk 


(2.5.2) 


where  y are  fixed  parameters  v*ith  a specified  structure.  ^ 

ljK 


E,tllk>  ■ »•  E<‘!Jk>  ■ »2 


E(eiJteiJ.>  ■ 0 V **' 


"'iJrW  - 0 °2"  ^k*  T*e 


E(eijrctks)  “ 0 p3*  J^k»  r*9- 


This  dispersion  matrix  of  (Y  fc)  can  be  exhibited  in  the  f- 


(2.1.3)  with  four  parameters  o , pj,  p^,  Pj.  A model  of  the  -vpe 


P 

■5, 


' * 


l i 


kirsi,  . . ^ Bj| 
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(2.5.2)  is  given  in  Rao  (1973,  p.  258). 


2.6  Mulcivsriats  Model 

A k-vari«te  linear  model  is  of  the  form 

(YjS...:^)  - X(B1:...,Bk)  + 

ECc^j)  - 0i])vi+‘**'Hjij)vp‘  (2.6.1) 

Denoting  Y - (Y| Y^),  ? - (8’ 6^).  7 - (ej e^) . the 

multivariate  model  may  be  written  as  a univariate  model 

Y - (I  © X)6  + t 

E(ce’)  ■ E (9.0V  ) (2.6.2) 

i-1  1 

where  are  (k*k)  matrices  of  variance  and  covariance  components. 

In  the  multivariate  regression  model  p • 1,  in  which  case 

E(EE')-0®V.  (2.6.3) 

We  may  specify  structures  for  e analogous  to  (2.2.1)  in  tii* 

univariate  case 


ci  " Vi^-^pV  1-1 k 

E(W]m,"a{j)'  E(*irV’°  T*B-  (2*6*4) 

For  special  choices  of  , we  obtain  multivariate  one,  two,... way 

mixed  models. 
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Models  of  the  type  (2.6.2)  have  been  considered  by  Krishnalah 
and  Lee  (1974).  They  discuss  methods  of  estimating  the  covariance 
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3.  ESTIMABILITY 


3.1  Unbiasedness 

Let  us  consider  Che  univariate  CM  model  (2,1.1)  with  the 
covariance  structure  (2.1.3) 


Y - XB  + €,  D(e)  - 0.V-+...+0  V 

11  p p 


and  find  the  conditions  under  which  linear  functions  f'0  can  be 
estimated  by  functions  of  Y subject  to  some  constraints.  The 
classes  of  estimators  considered  are  as  follows: 


(3.1.1) 


Q,  • {Y'AY,  A symmetric} 

Uf  - (g(Y):  E[g(Y>]  - f'0  V $eRm,  0cF) 
I - (g(Y)  - g(Y  + Xa)} 


(3.1.2) 

(3.1.2) 

(3.1.3) 


Theorem  3.1.1  provides  conditions  for  unbiased  estimabllity. 


Theorem  3.1.1.  Let  the  linear  model  be  as  in  (3.1.1).  Then: 


(1)  The  estimator  Y'AY  Is  unbiased  for  y ■ f'0  iff 
X'AX  - 0,  tr  AVf  - f , 1 - 1 

(11)  There  exists  an  unbiased  estimator  yeC)  Iff  feS(H), 

H - (htJ),  hAj  - trO^Vj  - PV^Vj) 

where  P is  the  projection  operator  onto  S(X). 

(Ill)  If  Y has  multivariate  normal  distribution,  then  Uff)  0 Is 
not  empty. 


(3.1.3) 


(3.1.6) 


3.1 


The  results  (1)  and  (11)  are  discussed  In  Seely  (1970),  Rao  (1970, 

1971)  and  Focke  and  Dowess  (1972)  and  (111)  In  Plncus  (1974). 

Note  li  Result  (11)  holds  If  In  (3.1.6)  we  choose 

hlJ  “ tr(Vl(I  " P)Vj)*  0.1.7) 

Note  2:  In  the  special  case  - 0 for  the  1-th 

Individual  parameter.  Is  unblasedly  estimable  Iff  t 0 where 

M - I - P. 


Lemma  3.1.1.  The  linear  space  P of  all  unblasedly  estimable 
linear  functions  of  9 i8 


T - (E  Y'AY:  A e spO^-PV^ V -PV  P)J  (3.1.8) 

6 P P 

where  sp(A. , . . . ,A  ) is  the  set  of  all  linear  combinations  of 
1 P 

A, , . . . ,A  • 

1 P 


Let  us  consider  the  multivariate  model  (2.6,2)  written  in  a vector 

form 


Y - (I  (C)X)0  + c 

E(ce')  - 0JOVJ+...+  Op0Vp  (3.1. 9) 

«**  matrix  variance-covariance  components. 

Lemma  3.1.2.  The  parametric  function  Y - £ f^  trC  0 is 
unblasedly  estimable  from  the  model  (3.1.9)  Iff  f'O  is  so  from 
the  univariate  model  (3.1.1). 


3.2 


Lemma  3.1.3.  The  class  r of  unblasedly  estimable  linear 


functions  of 


elements  of  0^,  1 ■ l,...,p  In  (3.1.9),  Is 


r ■ (y  ■ E trC^  are  such  that  Hb  ■ O^E  b^C^  ■ 0) 


where  H Is  as  defined  In  (3.1.6). 


(3.1.10) 


3 . 2 Invariance 

i 

An  estimator  Is  said  to  be  invariant  for  translation  of  the  parameter 
0 In  the  linear  model  (3.1.1)  If  it  belongs  to  the  class  (3.1.4). 

Theorem  3.2.1  provides  the  conditions  under  which  estimators  belonging 
to  the  class  ll^fl  I exist. 


Theorem  3.2.1.  Let  the  linear  model  be  as  in  (3.1.1).  Then: 


(i)  The  estimator  Y'AY  e l/f  A I iff 


AX  - 0,  tr  AVt  - fA,  i - 1 


(11)  There  exists  an  unbiased  estimator  in  class  Q_  A I iff 


i<V 


where 


(3.2.1) 


- (1'1  j ) » hAj  - tr (MV^MVj) , M - I - P. 


(3.2.2) 


(Hi)  Under  the  assumption  of  normality  of  Y,  the  result  (3.2.2) 
can  be  extended  to  the  class  I. 


Note:  In  (3.2.2),  we  can  choose 


hfJ  - tr (BB'V^BB'Vj) 

i 

where  B Is  any  choice  of  X , i.e.,  B is  a matrix  of  maximum 
rank  such  that  B'X  ■ 0. 


(3.2.3) 


•>  •» 


3.2.1.  The  linear  apace  of  all  invariantly  unblaaedly 


estimable  linear  functions  of  6 is 

Tj  - {E  Y’MAMY:  A e spO^-PV^, . . . ,Vp-PVpP)  }.  (3.2.4) 

0 

Lewsa  3.2.2.  If  f'0  is  invariantly  unbiasedly  estimable  from 
the  modal  (3.1.1)  then  so  is  y ■ I fj  trC  0t  from  the  model  (3.1.9). 

3.2.3.  All  invariantly  unbiasedly  estimable  linear  functions 
of  the  elements  of  0^...,  0,  in  the  model  (3.1.9)  belong  to  the 

•at 

p • (y  ■ I tr  Cj  0^:  Cj  are  such  that  Hj^b  ■ 0 Eb^C^  ■ 0).  (3.2.5) 

Note:  We  can  estimate  any  member  of  the  class  (3.2.5) 
by  functions  of  the  form 

Y - E tr (C j Y’AjY) 

where  Alt...,A  are  matrices  arising  in  invariant  quadratic  unbiased 
estimation  in  the  univariate  model  (3.1.1). 

3.3  Examples 

Consider  the  model  with  4 observations 

Yj  - $i  + Cj,  Y2  - &i  + e2*  Y3  * B2  + c3*  Y4  * 02  + C4 

2 

where  are  all  uncorrelsted  and  V(e^)  ■ V(c^)  " Oj  and 

V(c')  ■ V(e^)  • o*.  The  matrices  X,  Vj,  Vj  and  P the  projection 


operator  are  easily  seen  to  Vie 


X - 


V.  - 


0 

1 

0 

0 


0 

0 

0 

1 


O'1 

0 

0/ 


p - 


/i 

2 

1 

2 


0 


i i 0 


4 £ o 


o - — 

v n 


1°  0 i }/ 


The  matrices  H and  of  Theorems  3.1.1  and  3.2.1  are 


H -j 


3 

2 

1. 

2 


2 2 

Since  H is  of  full  rank,  and  are  individually  unblasedly 

estimable.  But  Is  of  rank  one  and  the  unit  vectors  do  not 

belong  to  the  space  SCH^)  and  therefore,  o^  and  are  not 

Individually  unblasedly  estimable  by  invariant  quadratic  forms. 

Consider  the  model  Y - XB  ♦ X$  + e where  S is  a fixed 

vector  parameter  and  $ Is  a vector  of  random  effects  such  that 

E<*>  - 0,  E(#4*>  - all  , F.(*e’)  - 0,  E(cc')  - o*I  . Let 

*■  m In 

2 

Y'AY  be  an  unbiased  estimate  of  o^.  Then  we  must  have 


X 'AX  - 0,  tr  AXX*  - 1,  tr  A • 0 


which  is  not  consistent  since  X'AX  - 0 tr  AXX'  - 0.  Hence 


unbiased  estimators  of  o \ do  not  exist. 


4.  MINIMUM  VARIANCE  UNBIASED  ESTIMATION-NORMAL  CASE 


4.0  Notations 

In  Section  3,  we  obtained  conditions  for  unbiased  estimabllity 
of  f'6  in  the  linear  model 

Y - XB  + e,  D(e)  - 6 V +...+6  V - V„ 

11  p p o 

restricting  the  class  of  estimators  to  quadratic  functions  of  Y. 

In  this  section  we  do  not  put  any  restriction  on  the  class  of 
estimators  but  assume  that 

y - Nn(xe,ve),  Mr”,  e c f 

i.e.,  n variate  normal,  and  Vg  is  p.d.  for  ® c F.  The 

condition  that  Vg  is  p.d.  is  assumed  to  simplify  presentation  of 

results,  and  is  satisfied  in  many  practical  situations. 

First,  we  derive  the  locally  minimum  variance  unbiased  estimator 

(IMVUE)  of  f'6  at  a chosen  point  ( B , 8 ) in  R^F . If  the 

o o 

estimator  is  independent  of  Bq,  0q  then  we  have  a uniformly 

minimum  variance  unbiased  estimator  (UMVUE) . Such  estimators  do 

not  exist  except  in  simple  cases.  In  the  general  case  we  suggest 

the  use  of  1MVUE  with  a suitable  choice  of  8 , ® based  on  previous 

o o 

experience  or  aprlorl  considerations.  We  also  Indicate  an  iterative 
method  which  starts  with  an  initial  value  ( 6^ , 0^ ) , gets  an  improved 
set  (B^,6j),  and  provides  in  the  limit  an  Iterated  MVE. 

LMVUE's  are  obtained  in  the  class  of  quadratic  estimators  by 
La  Motte  (1973)  under  the  assumption  of  normality  and  by  Rao  (1971) 
in  the  general  case.  Such  estimators  were  designated  by  Rao  as 


MIVQUE  (minimum  variance  quadratic  unbiased  estimator).  In  this 
section  we  show  that,  under  the  normality  assumption,  MIVQUF  is 
LMVUE  in  the  whole  class  of  unbiased  estimators. 

4.1  Locally  Minimum  Variance  Unhlased  Estimation 

Definition  4.1.1.  An  estimator  Y*  is  called  LMVUE  of  its 

expected  value  at  (0  ,6  ) c RmxF  iff 

o o 

V(Yj8o.eo>  < v<;|Bo,eo)  (4.1.1) 

for  all  y such  that 

E(Tt)'E(?)y(!,8)t  R^F.  (4.1.2) 


We  use  the  following  notations: 


v_  ■ e,v.+...+e  V 

0 11  p p 


Ai  ' ve  (vrpevip0)ve1-  pe  " Wv^xix’v;1 


k - (tr  aiVi) 


ka  a"  t(Y-X0)'A  (Y-X0) (Y-X  6)  'A  (Y-X6)]'. 

p,o  i P 


(4.1.3) 


Let  an  ^priori  value  of  (6,8).  Then  applying 

the  result  (3.1.6)  of  Theorem  3.1.1  we  find  that  f'8  is 
unblaaedly  estimable  iff 

f eS(r6  ). 


4.2 


Theorem  4.1.1  provides  an  explicit  expression  for  the  IMVUE. 


Theorem  4.1.1.  Let  f satisfy  the  condition  (4.1.4)  and  K0, 

^ . be  es  defined  in  (4.1.3).  Then  the  LMVUE  of  f’8  at  (0  ,8  ) 

8*8  o o 

is 

Y - X’k  - E X.(Y-X0  )*A  (Y-X0  ) (4.1.5) 

”o*  "o  o i o 

where  X Is  any  solution  of  ^ 1 • f. 

o 

Theorem  4*1«1  Is  established  by  shoving  that 
cov(g(Y),  y|80*«0)  - 0 

for  all  g(Y)  such  that  E[g(Y)|0,8j  » 0 V 0c  Rm,  8 c F,  and 
using  the  theorem  on  minimum  variance  estimation  given  in  Rao  (1973, 
p.  317). 

Hotels  For  any  X,  X'k.  is  LMVUE  of  its  expected  value 

“o’  ®o 

which  Is  a linear  functions  of  Thus  (4.1.5)  characterizes  all 

LMVUE* a of  linear  functions  of  8 at  0*  .9  )• 

o o 

Note  2;  The  variance  of  y as  defined  in  (4.1.5)  is 
V(y|0*0)  - 4(0-0o)'X*AVfl4X(0-0o)  + 2 tr  AV0AV0 
where  A ■ E XjA^  with  A^  computed  at  8q.  The  variance  at 

<V»„>  *• 

v(4|e0.e0)  - l n,i-i  f'*‘e  f 

o o 


(4.1.6) 


(4.1.7) 


4.3 


where  K is  any  g-inverse  of  K . 

v 0 

o o 


Note  3:  The  BLUE  of  X6  aC  6 is 
o 


X6  - P0  Y. 


Substituting  8 for  80  in  (4.1.5)  we  have 


(4.1.8) 


Ke  0 " ks  e 

Po  ®o’®o 


(4.1.10) 


Note  5:  Let  0 (i.e.,  each  component  of  0)  be  estimable  in 

which  case  K is  nonsingular  and  the  solution  of  (4.1.10)  is 
-1  ° 

91  " ^ a ’ Let  be  * solution  of  X8  - Pn  Y.  We  mav  use 

®o  bo,6o  1 0O 

®1*®1  LMVUE  of  6*8  as  initial  values  and  obtain  second  stage 

estimates  0j  and  8^  of  0 and  8 as  solutions  of 


;,.i'ks  - Y'(MV  M>+(U  V )<MV  M)\ 

o o o 

where  M ■ I - X(X'X)  X',  and  C-*"  is  Moore  Penrose  inverse  of  C 
(see  Rao  and  Mitra,  1972).  The  statistic  Yj  which  is  independent 
of  the  apriori  value  of  8 is  an  alternative  estimator  of  f'e  but 
it  may  not  be  unbiased  for  f'e. 


Note  4:  Theorem  4.1.1  can  be  stated  in  a different  form  as 

follows.  If  f'e  is  unblasedly  estimable  then  its  LMVUE  at  8 .0 

o’  c 

is  f'e  where  0 is  any  solution  of  the  consistent  equation 


(4.1.9) 


Kg  0 - kg  g , X6  - P*  Y. 
el  el,0l  01 


(4.1.11) 


4.4 


The  process  may  be  repeated  and  if  the  solutions  converge  thev 
satisfy  the  equations 


V - k6.e-  xe  - V- 


(4.1.12) 


The  solution  (0,9)  of  (4.1.12)  may  be  called  IMVUE  (iterated 
minimum  variance  unbiased  estimator)  of  (8,0).  The  exact  properties 
of  (8,0)  are  not  known. 


4.2  Invariant  Estimation 

Let  us  restrict  the  class  of  estimators  to  invariant  unbiased 
(IU)  estimators,  i.e.,  estimators  g(Y)  such  that 


g(Y  ♦ X8)  - g (Y)  V {* 
Efg(Y) |B.ei  - f'O 


(4.2.1) 


and  find  the  locally  minimum  variance  unbiased  estimator  (LMVIUE) , 


Let 


M - I - P,  P - X(X'X)“X' 

*\TI (®)  - (tr[(MVeM)+V1(MVflM)+V  ]) 


■ <tt'C(i-vvi-p;>%1v> 


hj(Y,0)  - (Y'(MVeM)+V1(W0M)+Y,...,Y,(MV0M)+V  (MVpM)+Y)' 


- (Y  -v;1  (i-p  ,)»,  o-ppvjh rv;1  (I-P  e)vp(i-r  • )v-,'y  i: 


(4.2.2) 


Theorem  4.2.1. 


(i)  f’0  is  invariantly  unbiasedly  estimable  iff 

4.5 


1 '1 


H 


pL 


n 

A 

J 


(4.2.3) 


f c *(^(0)) 

for  any  choice  of  0 such  that  V.  is  nonsingular , 


(ii)  The  LMV1UE  of  f’6  at  0 is 

o 


Y " ;,hT(Y,0o) 

where  X is  any  solution  of  [H^  (0  )Jx  - f. 


(4. 2.4) 


The  results  of  Theorem  4.2.1  are  obtained  by  transforming  the 
model  Y - XB  + e to  a model  involving  the  maximal  invariant  of  Y, 

t - B’Y  - R ' c - e 

* 

where  B ■ X^,  which  is  independent  of  6,  and  applying  Theorem  4,1,1. 

Note  1:  Theorem  4.2.1  can  be  stated  in  a different  form  as 
follows.  If  f'0  is  invariantly  unbiasedly  estimable,  then  its 
LMVIUE  at  0q  is  f’0  where  0 is  a solution  of 


(4.2.5) 


f«Ul(8o)j0  ' hI(Y'6o>  (4. 2. 6) 

where  Hur(0)  and  h^Y.G)  are  defined  in  (4.2,2). 

Note_2:  If  0 admits  invariant  unbiased  estimation,  then  as  in 
Note  5 following  Theorem  4.1.1  we  may  obtain  IMVIUE  of  (P,fl)  as  the  solution 
of 

xe  - peY  (4.2.7) 

fIlIII(0)  ]e  - hj  (Y,0). 


4.6 


5.  MINIMUM  NORM  QUADRATIC  ESTIMATION 
(MINQE-THEORY) 

5.0  MINQE-prlnclple 

In  Section  4 we  assumed  normal  distribution  for  the  random  vector 
Y in  the  linear  model  and  obtained  the  LMVUE  of  linear  functions  of 
variance  components  without  imposing  any  restriction  on  the  estimating 
function.  However,  we  found  that  the  estimators  were  all  quadratic. 

In  the  present  section  we  shall  not  make  any  distributional  assumptions 
but  confine  our  attention  to  the  class  of  quadratic  estimators  and  lay 
down  some  principles  for  deriving  optimum  estimators. 

Natural  estimator.  Consider  a general  random  effects  linear  model 
Y - X6  + U.*.+...+U  4>  - Xg  + U* 

11  p p 

E(^)  - 0,  E(^)  - 0^,  E(^j)  - 0,  i»«j.  (5.0.1) 

When  is  known,  a natural  estimator  of  0^  is 

®i  " Vi  + ri*  ri  " R<Gi>  (5.0.2) 

which  has  nice  properties  such  as  uniformly  minimum  variance  under 
mild  conditions  on  the  moments  of  Then  a natural  estimator  of 

* 

f'0,  a linear  function  of  0,  is 

f’0  - £(fi/r1)#i  - *'N*  (say).  (5.0.3) 

Given  a model  of  the  form  (5.0.2)  and  apriori  values  of  0^,  we 

can  by  a suitable  transformation  of  ^ 

♦l  ♦ ®i1/2Bi^i’  BiBi  * Gi* 


5.1 


MM 


write  the  model  in  the  form 


Y - XB  + U.*.+.. .+U  * 

11  p p 

1/2 

where  ■ U^a  ' 

E(*t)  - 0,  - etI  , E(^j)  - 0 (5.0.4) 

end  the  epriorl  values  of  0i  are  all  equal  to  unity.  We  assume  that 
such  a transformation  is  made  before  analysis  of  data.  In  terms  of  the 
model  (5.0.4)t  a natural  estimator  of  0^  is  ^^/r^. 

It  is  not  clear  how  the  concept  of  a natural  estimator  can  be 
extended  to  the  general  model 

Y - XB  + U<|>t  E(W’)  - eiFi+-**+epFp 

D(Y)  - eiUF1U'+...+epUFpU'  - ejVj+,.,+9  V (5.0.5) 

where  + may  not  have  the  structure  defined  in  (5.0.4).  However,  using 

-1/2 

prior  values  a^,...,ap  we  may  transform  $ to  F $ where 

1/2 

F “ aiFi+*  • ,+apFp*  and  11  to  WF  • (5.0.5)  represents  the  model 

after  such  a transformation  of  U and  + (0  and  Vj  remain  unchanged) 

we  may  formally  extend  the  concept  of  a natural  estimator,  although  we 
may  not  claim  any  optimal  properties  as  in  the  case  when  $ has  the 
structure  of  a random  effects  model. 

Definition  5.0.1.  Under  the  general  model  (5.0.5),  a natural 
estimator  of  y * f'O  is 

T • E w^'F^  - +'N+  (5.0.6) 


5.2 


where  N ■ I u.F  , p - (p,,...,u  )'  is  a solution  of  Hu  ■ f, 

11  i P 

H - (tr  F^). 

A more  general  definition  of  a natural  estimator  in  terms  of 
c when  the  model  is  Y ■ XfJ  + e without  specifying  any  structure 
for  e is  given  in  Section  5.4. 

MINQE-theory.  Consider  the  general  model  (5.0.5)  and  a quadratic 
estimator  y - Y'AY  of  f'6.  Now 


Y'AY 


• U'AU  U’AXW*. 

V '•X’AU  X'AX'V 


(5.0.7) 


while  the  natural  estimator  is  $'N$  as  defined  in  (5.0.6).  The 
difference  between  Y’AY  and  $'N$  is 


U’AU-N  U’AX.,*. 

V VAU  X’AX'V 


(5.0.8) 


The  minimum  norm  quadratic  estimator  (MINQE)  is  the  one  obtained  by 
minimizing  an  appropriately  chosen  norm  of  the  matrix  of  the  quadratic 
form  in  (5.0.8) 


D11 

D12 

U’AU-N 

U’AX 

°21 

D22 

X'AU 

X 'AX 

(5.0.9) 


We  shall  consider  mainly  two  kinds  of  norms,  one  a simple  Euclidean 
norm 

tr  DllDll  + 2 tr  Di2D21  + tr  D22D22  (5.0.10) 

and  another  a weighted  Euclidean  norm 


tr  D11WD11W  * 2 tr  D12KD?1W  + tr  D22KD22K 


(5.0.11) 


where  W and  K are  n.n.d.  matrices.  The  norm  (5.0.11)  gives 
different  weights  to  4>  and  (1  in  the  quadratic  form  (5.0.8). 

We  impose  other  restrictions  on  A (and  indicate  the  M1NQE  so 
obtained  by  adding  a symbol  in  brackets)  such  as  Y'AY 

(n)  is  unbiased:  MINQE(U) 

(b)  is  invariant  for  translation  in  8:  MINQF.(I) 

(c)  satisfies  both  (a)  and  (b) : MINQE(U,  I) 

(d)  is  unbiased  non-negative  definite:  MINQE (U,  0) 

(e)  is  invariant  non-negative  definite:  MINQE(I,  0),  etc. 

The  properties  of  the  estimator  strongly  depend  on  the  norm 
chosen  and  the  restrictions  imposed.  We  also  obtain  a series  of 
IMINQE's  (iterated  MINQE's),  by  repeatedly  solving  the  MINQE  equations 
using  the  solutions  at  any  stage  aa  prior  values  for  transforming  the 
model  as  indicated  below  equation  (5.0.5). 

5.1  MINQE (U,  I) 

We  consider  the  class  of  invariant  unbiased  quadratic  estimators, 
l.e.t  of  the  form  Y'AY  where  A belongs  to  the  class 

Cjj  - (A:  AX  - 0,  tr  AVt  - f , i-l,...,p)  (5.1.1) 

where  X and  are  as  defined  for  the  general  model  (5.0.5).  We 

use  the  following  notations  and  assumptions 

T - (V  + XX'J  > 0,  V - a,V,+...+a  V - UU' 
ci  all  p p 

PT  - x(x'T"1x)"x'T_1,  Rj  - (I  - PT) 

5.4 


Ml 


where  a Is  a prior  value  of  0. 

Theorem  5.1.1.  If  is  not  empty,  then  under  the  Euclidean 

norm  (5.0.10),  the  MINQE(U,  I)  of  f'6  is 

Y - I W*  Ai  " T_lRTViRTT~1 


where 

x - (x1§ 

. ,.,A  )’  is  any  solution  of 

(HUI(a)]X  - f 

where 

V°> 

is  the  matrix  (tr  A^V^). 

Proof.  Under  the  conditions  (5.1.1),  the  square  of  the 
Euclidean  norm  in  (5.0.10)  becomes 

IJU'AV-Wll2  - tr(U’AUU'AU)  - 2 tr  NU’AU  + tr  NN. 

But  N • I <>tr  NU’AU  - E p^f^  80  that  we  nee<*  minimire  only 

the  expression 

tr  AV  AV  ■ tr  ATAT  for  A e C,^T. 
o o UI 

It  is  easy  to  show  that  (5.1.5)  is  minimized  at  A*  such  that 

tr  DTAT  • 0V  I)  e C®  , 


(5.1.2) 


(5.1.3) 


(5.1.4) 


(5.1.5) 


(5.1.6) 


DeC®i  »>  D - R^,  E l^,  tr  E - 0 for  arbitrary  E.- 


Then 

T A*T  - E 


5.5 


which  gives  the  solution  (5.1.3).  The  equation  for  \ Is  obtained 

by  expressing  the  condition  of  unbiasedness.  Note  that  f (<x) J X ■ f is 

consistent  if  is  not  empty.  Also  the  solution  (5.1.2)  is  independent 

of  N. 

Note  Is  MINQE(U,  I)'s  are  additive. 

Note  2:  An  alternative  expression  for  y given  in  (5.1.3)  is 

Y “ Z AjY'A^,  - (MVaM)+V1  (MVM)+  (5.1.7) 

where  M • I - XX+. 

Note  3:  When  V is  nonsingular,  T can  be  replaced  by  V in 
a a 

Theorem  5.1.1.  Then 

' - £ Ai  ■ \\  viM  v;1.  (5.1.8) 

a a 

Note  4:  If  Y is  normally  distributed,  MINQE (U,  I)  is  LMVUE  of 
f'0  at  values  of  0 where  E 0 V is  proportional  to  V . (See 

1 1 Cl 

Theorem  4.1.1). 

Note  5:  If  in  (5.1.4)  we  use  the  weighted  Euclidean  norm  (5.0.11) 

| |U'AU-N| | 2 - tr(U'AU-N)  W(U’AU-N)W  (5.1.9) 

where  W is  p.d.,  the  solution  may  not  be  independent  of  N.  The 
expression  (5.1.9)  can  be  written  as 

tr  AGAG  - 2 tr  AH  + tr  NWN  (5.1.10) 

where  G ■ UWU'  and  H - UWNVU'.  If  G is  nonsingular,  then  the 
mipimua  of  (5.1.10)  is  attained  at 


5.6 


.>.-1 


a,  - e-‘«  rghrg*g 


(5.1.11) 


where  are  determined  from  the  equations 


tr  A^v^  ■*  f^,  i * 1».«.«P< 


Note  6:  It  is  seen  from  (5.1.2)  that  the  estimate  of  f'0 
can  be  written  in  the  form  f'0  where  0 is  a solution  of 


[HyT (a) ]0  - ho(Y,a) 


(5.1.12) 


where  the  i-th  element  of 


yY.a) 


is 


r»q>~ln  w O 'T-^V 

T 


VA.Y  « Y'T  RJ^KJT  -Y 


(5.1.13) 


and  Ha(U,l)  is  as  defined  in  (5.1.3).  If  each  component  of  0 

admits  invariant  unbiased  estimation  then  H (U,I)  is  non-singular 

a 

and  the  MINQE (U , I)  of  0 is 

0 - [HUI(o)]"1hI(Y,o)  (5.1.14) 

Note  7:  The  computation  of  MINQE (U,  I)  of  0 involves  the  use 
of  a a prior  value  o 0.  If  we  have  ..o  prior  information  on  0, 
there  are  two  possibilities.  We  may  take  a as  a vector  with  all  its 
elements  as  unity.  An  alternative  is  to  choose  some  a,  compute 
(5.1.14),  consider  it  (say  0^)  as  an  apriori  value  of  0 and  repeat 
the  computation  of  (5.1.14).  The  second  round  value,  say  0,,  is  an 
appropriate  estimate  of  0,  which  may  be  better  than  6^  if  the  initial 

a 

choice  a is  very  much  different  from  0^. 

A 

We  may  repeat  the  process  and  obtain  0^  choosing  0 ^ as  an 
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\ 


anriort  value  and  siVon.  The  limiting  value  which  satisfies  the 

/ 

equation 

[^^0)10  . hl(Y.e)  (5.1*15) 

in  the  IMINQE(U,  I),  the  iterated  MINQE (U , 1),  which  is  the  same 
es  IMVIUF.  defined  in  (4.2.7).  It  is  shown  in  Section  6 that  the 
equation  (5.1.15)  is  the  restricted  maximum  likelihood  (RML) 
equation  considered  by  Patterson  and  Thompson  (1975). 


5.2  MINQE (U) 

We  drop  invariance  and  consider  only  unbiasedness,  as  in 
problems  such  as  those  mentioned  by  Focke  and  Dewess  (1972)  where 
the  condition  for  invariance  does  not  hold.  In  such  problems  where 
Invariance  condition  is  not  used,  it  is  advisable  to  use  an  apriori 
value  of  0 and  change  Y to  Y - X6q  and  6 to  (F-8^) 

and  work  with  the  transformed  model  in  addition  to  the  transformation 
Indicated  in  (5.0.5).  The  class  of  unbiased  estimators  of  f'8  is 
defined  by 

Cj  - (A;  X’AX  - 0,  tr  AV^  - f^  i - 1 p}  (5.2.1) 

where  X and  are  as  in  the  general  model  (s  0.5). 


llieorc  i 5.2.1.  Let  - UU'  be  p.d.  If  is  not  empty 

then  the  MIIrQE(U)  under  Euclidean  norm  (5.0.10)  is 

Y - S At  - (Va  + XX’)"1  (V  - Py  V P^  ) (V^  + XX')'1 

a a 

where  X • (A^,.,.,^)'  is  any  solution  of 

[H(I(a)JX  - f 


(5.2.2) 


(5.2.3) 


where  Hy(o)  is  the  matrix 


(tr  A^). 


Proof.  Under  (5.0.10)  we  have  to  minimize 


I |U'AU-n| I2  + 2 | |u'AX| | 2 


which,  using  (5.2.1),  reduces  to 

tr  AV  AV  +2  tr  AV  AXX'_  tr  AV  A(v  + 2XX') 
a a a a a 

- tr  ATAT,  T - V + XX'. 

a 

The  expression  (5.2.4)  attains  a minimum  at  A^  iff 

tr  DTAT  - 0 V D c Q°. 

Observing  that  D c * F.  - PT?PPT  and  following  the 

arguments  of  Theorem  5.1.1,  the  expression  for  A^  is  obtained 
as  in  (5.2.2),  where  is  replaced  by  the  equivalent  expression 

PV  ‘ 
o 


(5.2.4) 

(5.2.5) 


(5.2.6) 


Note  1;  We  shall  consider  a few  alternatives  to  the  simple 

Euclidean  norm.  Focke  and  Dewess  (1972)  give  different  weights  to 

the  two  terms  in  (5.2.3)  as  in  (5.0.11).  Choosing  W ■ I and 
2 

K " r , (5.2.4)  becomes 

tr  AVaAVa+  2r2tr  AVaAXX\  (5.2.7) 

2 

The  constant  r determine  the  relative  weights  to  be  attached  to 
6 and  The  solution  obtained  by  minimizing  (5.2.7)  is  called 

r-MINQE(U)  which  is  the  same  as  (5.2.2)  with  X replaced  by  rX. 

* Note  2:  The  Iterated  estimates  of  0 and  MINQE(U)  of  0 are 
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solutions. of  the  equations 


x've‘x6  ■ x'v8It 

[Hy (8) ]0  - hyCY.e)  (5.2.8) 

where 

hyCY.e)  - (Y'A^Y Y'ApY)\  (5.2.9) 

Hy(0)  and  A^  are  as  defined  In  Theorem  (5.2.8).  The  solution 
of  (5.2.8)  is  represented  by  IMINQE(U). 

5.3  — MINQE(U) 

In  (5.2.7)  we  defined  r-MINQE(U)  which  uses  a weighted  Euclidean 
norm  to  provide  differential  wieghts  to  6 and  4 and  also  suggested 
a translation  in  Y using  a prior  value  of  8.  Actually  we  may 
consider  a transformation  which  changes 

-1  -1/2 

Y Y - X80,  8 -*■  r K ' 0 

2 

where  8q  and  r K correspond  to  aprlorl  mean  and  dispersion  of  8. 

Then  the  Euclidean  norm  of  (5.0.10)  becomes 

tr  A(Vq  + r2XKX’)  A(Va  + r2XKX')  (5.3.1) 

which  may  be  minimized. 

\ 

Let  ua  denote  the  optimal  solution  in  such  a case  by  A^  and 
define  Ag  - lim  Af  as  r •».  If  Ag  exists,  we  call  the 
correeponding  estimator  Y'AgY,  the  «-MINQE(U).  The  following  theorem 
duetto  Focke  and  Dewess  (1972)  establishes  the  existence  of  <»-MINQE(ll) . 
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Theorem  5.3.1.  Let  V be  the  set  of  linear  combinations  of 

(i)  «»-HlNQE(U)  exists  if  is  not  empty. 

(ii)  Aq  Is  the  unique  matrix  which  minimizes  tr  AVaAVa  in 

the  class  C - {A:  XKX'AV  M + MV  AXKX’  - V.-PV.P,  VA  e I/} 

a a 0 0 0 

(iii)  «®-MINQE(U)  is  invariant  with  respect  to  nonsingular  linear 
transformation  of  the  model  (5.0.5). 


Theorem  5.3.1  characterizes  «*-MINQE(U)  but  does  not  provide  the 
method  of  calculating  it.  Theorem  5.3.2  gives  the  formula  when  V - I, 
from  which  the  formula  for  general  V can  be  derived  by  a transformation 
of  the  model  in  view  of  statement  (ill)  of  Theorem  5.3.1. 


Theorem  5.3.2.  Let  G - (tr  MVjMV^),  B - (tr  MV^XKX') V ) and 

cf,  be  not  empty.  If  V - I,  the  «®-MINQE(U)  of  f'6  is  y ■ Y'A*Y 
U a 


(5.3.2) 


where 


A#  - (XKX')+V0M  + MV0(XKX’)+  + MVfeM 


vo-IaiV  vb-Ebivi 


(5.3.8) 


and  a ■ (a. , . . . ,a  ) ' , b - (b, , . . . ,b  ) ' satisfy  the  equations 
1 p x p 


Gb  + 2 Ba  ■ f,  Ga  - 0. 


(5.3.9) 


Theorem  5.3.3.  If  V is  p.d.  and  C..  is  not  empty,  the  «-MINQE(U) 

a u 

ia  obtained  by  replacing  M by  (MV  M)+  in  (5.3.8). 

a 

Proof.  We  consider  the  model  V~^2Y  - V-1^2X  + V_1^2U*  and 
a o a 

apply  the  result  of  Theorem  5.3.2.  On  slmplicatlon  the  formula  stated 
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in  Theorem  5.3.3  is  obtained. 


Note  1:  It  is  interesting  to  note  that  ~-MINQE(U)  is  the  same  if 

2 2 

instead  of  the  sequence  r K,  we  consider  (A  + r K)  for  any  A >_  0. 


Note  2;  --MINQE(U)  coincides  with  MINQE (U,  I)  if  it  exists. 


5. A MINQE  Without  Unbiasedr^ss 

Let  us  consider  the  linear  model 


Y - X3  + e,  E(cc')  - 0^1+. . .+epvp  - V0. 


Choosing  a prior  value  a of  0,  (5.A.1)  can  be  written 


1/2 

Y - XB  + V ' 

o * 


where  t.  * V c and  V **  a.V  +...+U  V . Using  the  definition 
* a a 1 p p 

(5.0.6)  with  e*  as  a natural  estimator  f'0  is 


(5.A.1) 


(5. A. 21 


r'tir  » f'ftl  V 1e  V ■ v-1^2 

c*we*  i i*'  **  vi*  a viva 


(5. A. 3) 


where  X ■ (X-.....X  )'  is  chosen  such  that  e'Ne.  is  unbiased 

1 p * n 

for  f'0,  i.e.,  X satisfies  the  equation  [H(o)]X  ■ f where 
H(«)  - (tr  Vt,»  ) - (tr  V^VjV’Vj). 


(5. A. A) 


It  is  seen  that  (5. A. 3)  is  LMVUE  of  0 at  0 ■ a when  e is 
normally  distributed.  In  most  of  the  applications  the  natural 
estimator  is  independent  of  a,  which  is  an  ideal  situation. 
Thfe  MINQE  of  f'0  is  Y'AY  where  A Is  chosen  to  minimize 


5.12 


(5.4.5) 


v1  "av 1 ' 2-n 

a a 

1 >2 

X'AV  ' 
a 

1/2 

V ‘ax 

a 

X ' AX 

In  Sections  5.1  - 5.1  we  imposed  the  condition  of  unbiaseines; 
on  Y AY.  We  withdraw  this  condition  but  consider  some  altorn.it'  *• 
restrictions  on  the  symmetric  matrix  A as  defined  bv  the  follow'd"' 

classes. 


C - {A} 

(5.4.6) 

Cpu  - {A:  X'AX  - 0) 

(5.4.7) 

Cj  - (A:  AX  - 0) 

(5.4.8) 

It  is  seen  that  when  ArCpjI,  the  bias  in  the  estimator  Y 'AY  is 
Independent  of  the  location  parameter  R,  and  is  thus  partial lv 
unbiased  (PU) . The  MINOR ' s obtained  subject  to  the  restrict  tor.'- 
(5.4.6)-(5.4.8)  are  represented  by  MINQE,  MINQE(PU),  MINOE(T> 
respectively.  The  concept  of  MINQE(I)  is  due  to  Poduri  Rao  and 
Cbaubey  (1978). 


Theorem  5.4,1.  Consider  the  model  (5.0.5)  and  let 

Va  ” °tlVl'K,-+Vp  1 V ch<~.  1 c W " S AtV  where 

X " (Xj»...tX  )'  satisfies  the  equation  (H(a)]X  *«  f,  where 
^ («)  " ^ i^cx  Xr  ^ • Then  under  the  Euclidean  norm  in 

(5.4.5),  the  optimal  matrix  A*  providing  MINQE' s are  as  follows. 

(i)  MINQE:  A - (V  + XX')-1W(V  + XX')*1 

* a a 

(it)  MINQE  (Pit)  : A.  - (V  i XX')~l(W  - P WP  ) (V  + XX')-' 

« ana 

P - X(X'V  x)"x'v"1 
o an 


(5.4.9) 


> 5.4.10) 


5.13 


(ill)  MINQE(T) : A - (MV  M)  W(MV  M)  , 


V~\l  - P)W(I  - P iv’1 
a a a a 


where  M - I - X(X'X)“X' 


Proof.  Under  Euclidean  norm,  the  square  of  (5.4.5)  is 


1/21/72  7 

tr(V  ' AV  -N)  + 2 tr (X'AV  AX)  + tr(X'AX)  , 

a a a 


Without  any  restriction  on  A,  the  minimum  of  (5.4.12)  is  attained 


at  Aa  iff 


tr(V1/2A^V1/2-N)V1/;RV1 12  + 2 tr(X'A*V  BX)  + tr(X'A.XX’BX)  *=  0 


a a 


for  any  symmetric  matrix  R.  Then  A satisfies  the  equation 


V1/2(V1/2A.V1/2-N)V1/2  + XX'A.V  + V A. XX ' + XX'AJCX'  = 0 
a a * a ex  * a a * * 


or 


l/2„fl/2 


(V  + XX ' )A . (V  + XX')  - V ' NV ' *=  Z A V = w 


a a 


i i 


A.  = (V  + XX')  W(V 
* a a 


')_1W(V  + XX')  1 


which  is  the  matrix  gi.en  m (5. a. 9). 


If  A is  subject  to  the  restriction  X'AX  * 0,  then  (5.4.13) 


must  hold  when  B is  replaced  by  B - P'BP  where  P is 

a a a 


defined  in  (5.4.10) . Then  arguing  as  above  and  noting  that  TV  *>  V P 

an  a i 


the  equation  for  is 

(V  + XX' ) A . (V  + XX')  * F.  X . (V , - P V P') 
'a  * a il  ala 


(5.4.11) 


(5.4.12) 


(5.4.13) 
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{ 

\ 

or  A.  - (V  + XX')_1(W  - P WP')(V  + XX')-1 

, * a a a a 

which  Is  the  matrix  given  in  (5.4.10). 

If  A is  subject  to  the  condition  AX  ■ 0,  then  (5.4.13) 
must  hold  when  B is  replaced  by  MBM  where  M « I - P.  Then  A^ 
satisfies  the  equation 

(MV  M)A  (MV  M)  - MWM 
a * a 

or  A^  - (MVoM)+W(MV  M)+ 

« V-1(I  - P )W(I  - P'JV-1 
a cr  ' a a 

which  is  the  matrix  given  in  (5.4.11). 

Note  1:  MINQE  in  (5.4.9)  and  MINQE(I)  in  (5.4.11)  are 
automatically  non-negative,  while  MINQE (PU)  may  not  be. 

Note  2;  The  MINOE(I)  of  f’0  given  in  (5.4.11)  can  be  written 

a A 

as  f’9  where  0 is  a solution  of 

(5.4.14) 


(5.4.15) 


[H(a)]0  = h^Y.a) 

where  H(a)  is  as  defined  in  (5.4.4)  and  the  i-th  element  of 
hi(Y,a)  is 

Y,V~1(I  - P )V  (I  - P')V-1Y. 
w ct  i a a 

The  equation  (5.4.14)  is  consistent.  If  0 is  identifiable  then 
H(0)  is  non-singular,  in  which  case  0 * fH(o)]-1h  (Y,I).  This 

(t 

form  of  the  solution  pnab.1'''~  us  to  obtain  IMINQE(I),  i.e.,  iterated 
0 
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MINQE(I),  by  writing  9^  - [H (a ) ] *ha(Y,I)  and  obtaining  a second 
stage  estimate  6 with  a replaced  by  9^.  The  limiting  solution, 
if  the  process  converges,  satisfies  the  equation 

[H(»)]e  = h^Y.Q)  (5.4.16) 

which  is  shown  to  be  the  maximum  likelihood  equation  in  Section  6. 


5.5  MINQE(D)-Non-negative  Definite  Estimator 

In  the  general  variance  components  model,  we  admitted  the 
pos88iblllty  of  some  of  the  parameters  being  negative.  But  there 
are  cases  such  as  the  random  effects  model  where  the  variance  components 
are  non-negative  and  it  may  be  desirable  to  have  non-negative  estimators 
for  them.  The  estimators  considered  so  far  except  those  in  Section  5.4 
can  assume  negative  values  although  the  parametric  function  is  non- 
negative.  In  this  section  we  explore  the  possibility  of  obtaining 
unbiased  quadratic  estimators  y * Y'AY  with  A 0 of  parametric 
functions  f'0  which  are  non-negative  in  6eF  for  a general  model. 

A MINQE  in  this  class  is  denoted  by  MINQE(U,  D) , where  D stands  for 
non-negative  definiteness  cf  the  quadratic  estimator. 

The  following  lemma  characterises  the  nature  of  the  matrix  A if 
Y has  to  be  unbiased  and  non-negative. 

Lemma  5.5.1.  A non-negative  and  unbiased  quadratic  estimator 
Y'AY  satisfies  the  invariance  condition,  i.e.,  AX  - 0. 

Proof.  Unbiasedness  r^X' AX  - Os^AX  » 0 since  A >_  0 . 

* In  view  of  Lemma  5.5.1  we  need  only  consider  the  class  of  matrices 
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i 


CjD  - {A:  A > 0,  AX  - 0,  tr  AVf  - f±,  i-l,...,p). 


(5.5.1) 


Further , because  of  invariance  we  can  work  with  a transformed  model 


t - Z'Y  - e 


E(t)  - o,  E(tt')  - e1B1+...+epBp 


(5.5.2) 


where  Z - X (with  full  rank  say  s)  and  Bj^  * Z'V^Z,  i - l,...,p. 
We  need  consider  quadratic  estimators  y - t'Ct  where  C belongs  to 
the  class 


ciro  " {C:  c - °»  tr  CBi  * V* 


(5.5.3) 


Lemma  5.5.2.  C*D  is  not  empty  iff 


f c convex  span  (q(b):  beRn) 


(5.5.4) 


where  q(b)  - (b'MVjMb b'MV  Mb)', 


Note;  In  terms  of  the  model  (5.5.2),  the  condition  (5.5.4)  in 


f c convex  span  (q(b),  bcR8} 


(5.5.5) 


where  q(b)  ■ (b 'B,b, . . . ,b 'B  b) . 

1 P 


The  conditions  (5.5.4)  and  (5,5.5)  Asa  rather  complicated,  but 
simple  results  can  be  obtained  if  we  assume  V to  be  n.n.d. 


Theorem  5.5.1.  Let  >_  0,  i - l,...,p,  V ■ E and  ^ V-V^ 

and  be  as  defined  in  (5.5.2).  There  exists  an  n.n.d.  auadr.itir 


unbiased  estimator  of  6 . iff 
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1 1 
I $ 


S(B ^)<h  s (*  1^)4*  SCON ^ <£  5(MV(1)M) 


I 


4=>  5(MV(j)M)^>  S(MVM)4->R(MV(1)M)  < R(MVM)  . 


Note  X:  The  condition  (5.5.6)  can  also  be  expressed  as 

(I  - GJVjd  - G)  * 0 (5.5.7) 


(5.5.6) 


where  G is  the  projection  operator  onto  the  space  Renerated  by  the 
columns  of  the  compound  matrix 


(X: 


r 


(5.5.8) 


Note  2:  If  S^JDSfM),  then  S(MVjM)  P 5(MVtM)  for  all  1, 
In  which  case, application  of  Theorem  5.5.1  shows  that  at  most  Is 

non-negatlvely  estimable. 


Not*  3:  If  S(V])  D S(M)  and  S(V2)  3 S(M),  then  none  of  the 
single  components  are  non-negatlvely  estimable. 

Note  4:  (LaMotte.  1973.)  If  V(J)  > 0,  then  6,  Is  not  non- 
negatlvely  estimable.  Further,  If  V1  > 0,  then  , i^j  is  not 

non-negatlvely  estimable. 

However,  let  us  assume  that  C*D  is  not  empty  for  a given  f and 
estimate  f’e  by  MINQF  principle.  For  this  purpose  we  have  to  minimize 

||A||~  - tr  AVAV  when  A cC^  . (5.5.9) 


* 
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* 

m . 


This  appears  to  be  a difficult  problem  in  the  general  case.  Of  course, 
if  MINQE(U,  I)  turns  out  to  be  a non-negative  estimator  in  any  given 
situation  it  is  automatically  MINQF.(ll,  D) . It  may  also  be  noted  that  if 
sp{MVjM, . . . ,MVpM}  is  a quadratic  subspace  with  respect  to  (MVM)+t 
then  the  MINQE(U,  I)  of  f'6  is  n.n.d.  iff  is  not  empty. 

Since  C^D  is  a convex  set,  we  proceed  as  follows  to  solve 


the  problem  (5.5.9).  The  minimum  Is  attained  at  Aa  iff 

tr  BVAV  > tr  A*VA*V  y BeC^  (5.5.10) 

or  writing  B ■ A*  + D,  the  condition  (5.5.12)  becomes 

tr  DVA*V  > 0 V DeP  (5.5.11) 

V - {D:  DX  - 0,  A*  + D >_  0,  tr  - 0,  i - l,...,p}.  (5.5.12) 


A general  solution  for  (5.5.11)  cannot  be  explicitly  written  down,  but 
the  formula  will  be  useful  in  examining  whether  any  guessed  solution 
for  A^  provides  a MINQE(ll,  D) . We  shall  consider  some  special  cases. 

Theorem  5.5.2.  Let  >_  0,  i ■ l,...,p,  and  0^  be  estimable,  i.e., 
the  condition  (5.5.0)  (s  satisfied.  Then  the  MTNQE(I’.  0)  of  0^  is 

“ pXA^T  Y’AiY*  Aj  " (<T  ' C^Vj(I  “ G)l+  (5.5.13) 

where  G is  the  projection  operator  onto  the  space  generated  by  the 
columns  of  (X,  Vj,...,Vj_j,  V Vp) . 

An  alternative  approach  to  the  problem  (5.5.9)  based  on  standard 
methods  of  convex  programming  1=  provided  by  Pu  kelshelm  (1077). 
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We  define  the  functional 


g(B)  - min  (||a||  -*A,B>) 

AcCL 


(5.5.14) 


2 

where  | |a|  | ■ tr  AVAV  and  <A,B>  - tr  AVBV  with  V > 0, 


call  the  problem 


sup  g (B) 

B>0 

as  the  dual  optimization  problem. 

Lemma  5.5.3.  Let  A^eC^  and  B*  >_  0 be  such  that 

IIaJI2  - g(Bj 


(5.5.15) 


Then: 


(5.5.16) 


(i)  A^  and  B^  are  optimal  solutions  of  (5.5.9)  and  (5.5.15) 
(ii)<A*,B*>  -0. 

Note:  g(B)  is  bounded  above  since 

|!aJ|  >_*(B)  for  all  B. 


For  obtaining  a satisfactory  solution  to  the  problem  (5.5.9)  we 
need  an  explicit  expression  for  g(B).  We  obtain  this  in  terras  of  A 
where  Y AY  is  the  MTNQF(li,  I)  of  f'6.  Let  us  note  that  any  matrix 
B(*B')  can  be  decomposed  in  terms  of  symmetric  matrices 


B - flU  + (B  - B°) 


such  that  B°e(’lj  and  < Bl  , o.  The  matrix  B°  is  simplv  the 


(5.5.17) 


(5.5.18) 


M 
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projection  of  B onto  the  subspace  in  the  space  of  symmetric 

matrices  with  inner  product  < *,*>  as  defined  in  (5.5.15).  We  note 
that  by  construction  A is  such  that 

„ o 

<A,  B > » 0 for  any  Riven  B. 


Theorem  5.5.3.  I.et  Y'AY  be  MINQE(ll,  I)  of  f'9  and  be 

not  empty.  Then: 


(5.5.19) 


(i)  g(B)  - ||A|!2  -<A,B>  - J | | R°  | | 2 


(ii)  B^  is  optimal  (i.e.,  maximizes  g(B)]  iff 


A + ; R®  > 0,  < A + B*.  B*>  - 0 


(iii)  A*  - A + | B* 

is  a solution  to  (5.5.9),  i.e.,  provides  MINQE(U,  D)  of 
f'0  and  < A . ,B.  > = 0. 

a ft 


(5.5.20) 


(5.5.21) 


(5.5.22) 


The  results  of  Theorem  5.5.3  are  still  complicated.  A sufficient 
condition  for  optimality  of  A*  is  given  in  Theorem  5.5.4. 


Theorem  5.5.4.  If  there  exists  a B > 0 such  that  and 


A*  " A 


T«"i. 

I I B°  2 


then  A^  is  an  optimal  solution  of  (5.5.11) 


Proof . Check  optimality  t>v  the  condition  (5.5.1b)  noting  that 


R «•  -2  — , B. 

* I I y.0  | | 2 
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6.  MAXIMUM  LIKELIHOOD  ESTIMATION 


6.1  The  General  Model 

We  consider  the  general  GM  model 

Y = XB  + e 

E(ee')  ■ 0,V  +...+0  V -V  (6.1.1) 

11  P P 0 

and  discuss  the  maximum  likelihood  estimation  of  0 under  the 
assumption 

Y - N (XB ,V  ) , BeRm,  0cF.  (6.1.2) 

n 0 

We  assume  that  V is  p.d.  for  V 0eF. 

Harvllle  (1977)  has  given  a review  of  the  ML  estimation  of  0 
describing  the  contributions  made  by  Anderson  (1973),  Hartley  and 
Rao  (1967),  Henderson  (1977),  Patterson  and  Thompson  (1975),  Miller 
(1977,  1979)  and  others.  We  discuss  these  methods  and  make  some 
additional  comments. 

The  log  likelihood  of  the  unknown  parameters  (6,0)  is 
proportional  to 

*(B,0,Y)  = - log j Vg | - (Y-XB)'  V“1(Y-XB).  16.1.3) 

A A 

The  proper  ML  estimator  of  (6,0)  is  a value  (6,0)  such  that 


i(M,Y)  - sup  £(8,0,Y).  (6.1.4) 

8,0eF 

Such  an  estimator  does  not  exist  in  the  important  case  considered 
by  Focke  and  Dewess  (1972).  In  the  simple  version  of  their  problem 
there  are  two  random  variables 


(6.1.5) 


Yl  " ''  + Gl*  Y2  " M + C2 


E(eJ)  - o*.  E(e2)  = a2,  E (e^)  - 0. 


The  likelihood  based  on  Y, 


and  Is 


(YyyV 

- log  0l  - log  o2  - -y- 


2o 


(y2-w)‘ 

2a2 


(6.1.6) 


which  can  be  made  arbitrarily  large  by  choosing  u « Y^^  and  letting 

-►  0,  so  that  no  proper  MLE  exists.  The  ML  equations  obtained  equating 
the  derivatives  of  (6.1.6)  to  zero  are 

o2  - (Yj-y)2,  cr2  « (Y2~y)2,  y(-y  + -y)  ■ -y  + 


_2 

°1  °2  °I 


(6.1.7) 


which  imply  o^  » o2»  Thus  the  ML  approach  fails  to  provide 

acceptable  estimators.  However,  in  the  example  (6.1.5),  all  the 

2 2 

parameters  are  identifiable  and  MINQE(U)  of  o^  and  cy  exist. 

2 2 

A similar  problem  arises  in  estimating  and  o 2 in  the  model 

Y « XP  + Xy  + e where  F(yy')  ■ o2I  , F.(ee')  ■>  o2l  and 

1 m (in 

F(fy’)  " 0. 

It  is  well  known  that  hi.  estimators  of  variance  components  ore 
heavily  biased  in  general  and  in  some  situations  considered  by  Neyman 
and  Scott  (1948),  thev  are  not  even  consistent.  In  such  cases,  the 
use  of  ML  estimators  for  drawing  inferences  on  individual  parameters  may 
lead  to  gross  errors,  unless  the  exact  distribution  of  the  Ml,  estimators 
is  known.  These  drawbacks  and  the  computational  difficulties  involved 
in  obtaining  the  Ml,  estimators  make  the  ML  method  less  attractive  for 


practical  applications. 


6.2  Maximum  Likelihood  Equations 


For  OeF  such  that  V > 0 (i.e.,  p.d.),  the  likelihood  of 

0 

(B,0)  la 

1(0, 0,Y)  - - log|V0|  - (Y-X0)  ’V^CY-XB)  . 

Taking  derivatives  of  (6.2.1)  w.r.t.  to  0 and  0^  and  equating 
them  to  zero  we  get  the  ML  equations 

X'V^XB  - X'V^Y 
0 0 

tr  V;1V1  - (Y-X0)’V“1V1V"1(Y-X0) 

1 - 1 , . . . ,p. 

Substituting  for  0 In  (6.2.3)  from  (6.2.2),  the  equations  become 

X6  * p0y,  P0  - X(X'V01X)"X'V01 

fH(0)]0  - hj(Y,e> 

where  H(0)  » (tr  V0*V^V0*Vj)  is  the  matrix  defined  in  (5.4.4)  and 
the  i-th  element  of  hT(Y.9)  Is 

yd  - %yv;\v;la  - p6)y 

which  la  the  same  as  the  expression  defined  in  (5.4.15). 

We  make  a few  comments  on  the  equations  (6.2.4)  and  (6.2.5). 

(1)  The  MI.  equation  (6.2.5)  is  the  same  as  that  for  TMTNQF.(I) 
given  in  (5.4 .15) . 

(11)  The  original  likelihood  equation  (6.2.3)  Is  unbiased  while 
the  equatir..,  , which  provides  a direct  estlncie  of 


(6.2.1) 

(6.2.2) 

(6.2.3) 

(6.2.4) 

(6.2.5) 

(6.2.6) 
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is  not  so  in  the  sens? 


Elh  (Y,6)]  rH(0)]0.  (6.2.6) 

An  alternative  to  the  equation  (6.2.5)  is  the  one  obtained 
by  equating  h^(Y,0)  to  its  expectation 

hjOr,e)  - Ffh^Y.Q)]  - [HUI(9)]0  (6.2.7) 

which  is  the  restricted  ML  (RML)  equation  suggested  by 
Patterson  and  Thompson  (1975) . 

(iii)  There  may  be  no  solution  to  (6.2.5)  in  the  admissible  set 

F to  which  0 belongs.  This  may  happen  when  the  supremum 
of  the  likelihood  is  attained  at  a boundary  point  of  F. 

(iv)  It  is  interesting  to  note  that  the  ML  estimate  of  0 is 
invariant  for  translation  of  Y by  Xo  for  any  n,  i.e., 
the  MLE  is  a function  of  the  maximal  invariant  B'Y  of  Y 
where  B *■ 

Suppose  0 in  the  model  (6.1.1)  is  identifiable  ov 
the  basis  of  distribution  of  Y in  the  sense: 

0. V.+...+0  V - o 'V.+...+0'V  «=,>  0,  - 0’  = 0 for  o.i  i. 

11  p p I 1 p p i i 

1. e. . V ■ r.  1 • f • vr  ’ • Independent.  Put  it  mr»y  ’ 

as  in  the  case  of  the  example  of  Focke  and  ftewess  (1972), 
that  0 is  no  longer  identifiable  when  we  consider  only 
the  distribution  of  B'Y,  the  maximal  invariant  of  Y. 

Such  a situation  arises  when  R’V^R  are  linearly  dependent 
while  are  not.  In  such  cases  the  ML  method  is  not 

applicable  while  MTNOF(U)  developed  in  Section  5.2  ran  be 
used.  Thus,  the  invariance  property  of  MLF  limit'-  t be 
scope  of  application  of  the  ML  method. 


(v)  Computational  algorithms:  The  equation  (6.2.5)  for  the 
estimation  of  9 is,  in  general,  very  complicated  and 
no  closed  form  solution  is  possible.  One  has  to  adopt 
iterative  procedures.  Harville  (1977)  has  reviewed  some 
of  the  existing  methods. 

(a)  If  0^  is  the  k-th  approximation  to  the  solution  of 
(6.2.5),  then  the  (k+l)-th  approximation  is 


Vi  ■ 


(6.2.8) 


as  suggested  for  IMINQE(I),  provided  6 is  identifiable, 
on  the  basis -of -the  maximal  invariant  t>6 — V.  Otherwise, 
the  H matrix  in  (6.2.5)  is  not  invertible.  Iterative 
procedure  of  the  type  (6.2.8)  is  mentioned  by  Anderson 
(1973),  Harville  (1969),  LaMotte  (1973)  and  Rao  (1972) 
in  different  contexts.  However,  it  is  not  known  whether 
the  procedure  (6.2.8)  converges  and  provides  a solution 
at  which  supremum  of  the  likelihood  is  attained . 

(b)  Hartley  and  Rao  (1967),  Henderson  (1977)  and  Harville 
(1977)  proposed  algorithms  suitable  for  the  special  case 
when  one  of  the  is  an  identity  matrix  (or  at  least 

non-singular) . An  extension  of  their  method  for  the  general 
case  is  to  obtain  the  (k+l)-th  approximation  of  the  i-th 
component  of  0 as 


0 = 0 
i,k+l  ik 


y’u-p;  rvsVv^a-p,  )y 

ek  ek  1 ek  ek 


tr  VjHr, 


(6.2.9) 


i ■ 1 > • • • » P • 


In  the  special  case  when  are  non-negative  definite 


b 


and  the  initial  0^  are  chosen  as  non-negative,  the 
successive  approximations  of  0^  using  the  algorithm 
(6.2.9)  stay  non-negative.  This  may  be  a "good  property" 
of  the  algorithm,  but  it  is  not  clear  what  happens  when 
the  likelihood  equation  (6.2.5)  does  not  have  a solution 
in  the  admissible  region. 

(c)  Hemmerle  and  Hartley  (1973)  and  Goodnight  and  Hemmerle  (1978) 
developed  the  method  of  W transformation  for  solving  the 

ML  equations.  Miller  (1979)  has  given  a different  approach. 

Possibilities  of  using  the  variable-metric  algorithms  of 

Davidson-Fletcher-Powell  described  by  Powell  (1970)  are 

mentioned  by  Harville  (1977).  As  it  stands,  further 

research  is  necessary  for  finding  a satisfactory  method 

of  solving  the  equation  (6.2.5)  and  ensuring  that  the 

solution  provides  a maximum  of  the  likelihood. 

I 

6.3  Restricted  Maximum  Likelihood  Equation 

As  observed  earlier  the  ML  equation  (6.2.5)  is  not  unbiased,  i.e., 

E l h j (Y , 0 ] i [H(0)]0.  (6.3.1) 

If  we  replace  the  equation  (6.2.5)  by 

hI(Y,0)'-  E[hj (Y,0) ] 

- [Hux(0)]0  (6.3.2) 

we  obtain  the  IMINQE(U,T)  defined  in  (5.1.14),  which  is  the  same  as 
IMVIUE  defined  in  (4.2.7). 

The  equation  (6.3.2)  Is  obtained  by  Patterson  and  Thompson  (1975) 

0 

by  maximizing  the  likelihood  of  0 based  on  T'Y,  where  T is  any 
choice  of  X , which  is  the  maximal  invariant  of  Y.  Now 


6.6 


(6.3.3) 


£ (0 ,T ' Y)  = - lop  | T ' V T I - Y 'T  (T ' V T)  T’Y. 

0 o 


Differentiating  (6.3.3)  w.r.t.  0^  we  obtain  the  RML  (restricted 


ML)  equation 


tr(T(T'V0T)-1T'V.)  - Y 'T (T 'V  T)~^T 'V  .T (T 'V  T) _1T 'Y 
0 1 0 10 


(6.3.4) 


i **  1 | • • • y P I 


Using  the  identity  (Rao,  1973*  p.  77) 


T (T'V  T) _1T ' - V~T  - V~1X(X,V"1X)X’V"1 
0 0 0 0 0 


■Ve'(I  -V 


(6.3.5) 


the  equation  (6.3.4)  becomes 


^eVw  - Vv-^i-p^v^i-p^v^ 


(6.3.6) 


i ■ 1 > . . . «p 


which  la  Independent  of  the  choice  of  T ■ X used  in  the  construction 
of  the  maximal  invariant  of  Y.  It  is  easy  to  see  that  (6.3.6)  can  be 

written  as 

[Hux(0)]0  - hi(Y,0)  (6.3.7) 

which  is  the  equation  (6.3.2). 


(i)  Both  ML  and  RML  estimates  depend  on  the  maximal  invariant 
T'Y  of  Y.  Both  the  methods  are  not  applicable  when  0 
is  not  identifiable  on  the  basis  of  t'Y. 

(ii)  The  bias  in  RMI.E  may  not  be  as  heavy  as  in  MLE  and  may  be  more 
useful  as  point  estimators. 

(ill)  The  solution  of  (6.3.7)  may  not  lie  in  the  admissible  set 
of  0 as  in  the  case  of  the  ML  equation. 


. 6.7 


* . 

. ■< 


p 


R 


■ i 


— .I 


(iv)  If  0^  is  the  k-th  approximation,  then  the  (k+l)-th 
approximation  cun  he  obtained  as 


Vi  • *"6  «'.t>rV<Y,D. 

k k 


(6.3.8) 


It  is  not  known  whether  the  process  converges  and  yields 
a solution  which  maximizes  the  marginal  likelihood. 

(v)  Another  algorithm  for  RM1.F.  similar  to  (6.2.9)  is  to  compute 
the  (k+l)-th  approximation  to  the  i-th  component  of  0 as 


9i,k+l 


i.k 


Y'a-p5  )»jV*j’(i-pj )» 

k k k k 


(6.3.9) 


It  is  seen  that  both  Ml.  and  KML  estimators  can  be  obtained  as 
iterated  MINQE's,  MI.E  being  IMINQE(I)  defined  in  (5.4.16)  and  RMLE 
being  IMINQE(U,  I)  defined  in  (5.1.14).  There  are  other  iterated 
MINQE's  which  can  be  used  in  cases  where  Ml.  and  RM1,  methods  are  not 

applicable. 

It  has  been  remarked  bv  various  authors  that  MTNQF  Involves 
heavy  computations,  requiring  the  Inversion  of  large  matrices.  This 
argument  is  put  forward  against  the  use  of  MINQE.  These  authors 
overlook  the  fact  that  inversion  of  large  matrices  depend  on  the 
Inversion  of  smaller  order  matrices  in  special  cases.  For  instance, 
if  Vg  is  of  the  form  (T  + UPU’),  then  it  is  well  known  that 

V(,1  - 1 - lT(tTU  + trW  (6.3.10) 

which  can  be  used  to  compute  If  the  matrix  (IT  *1’  + P * ) is 

compara tivel v of  a smaller  order  than  V . it  mav  be  noted  that  the 

computational  comp  1 o'  l\  i * 1.  * same  eider  fei  UNQF  and  v!.|-. , RMl.E. 


Whan  ve  have  good  prior  Information  MINQE's  should  be  better  than 

•••  • 

MLE’e  or  RMLE’s. 
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